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The centrality dependence of y/s NN = 200 GeV d+Au J/ip data, measured in 12 rapidity bins 
that span —2.2 < y < 2.4, has been fitted using a model containing an effective breakup cross section 
combined with EPS09 shadowing. The centrality dependence of the shadowing contribution was 
allowed to vary nonlinear ly, employing a variety of assumptions, in an effort to explore the limits of 
what can be determined from the data. It is found that the onset of shadowing is a highly nonlinear 
function of impact parameter. The impact parameter dependencies of the effective breakup cross 
section and the shadowing parameterization are sufficiently distinct to be determined separately. 



PACS numbers: 25.75.Dw 

The modification of the gluon distributions in nuclear 
targets in high energy collisions, referred to here as gluon 
shadowing, is inherently interesting because of what it 
can teach us about the behavior of gluons at low Bjorken 
momentum fraction, x, where the gluon densities are high 
and saturation effects are expected to become impor- 
tant pQ. In addition, the modification of parton distribu- 
tions in nuclei determines the initial conditions in a high 
energy nuclear collision. The initial conditions must be 
sumcently well understood before final-state hot matter 
effects can be isolated. 

Parameterizations of the dependence of nuclear- 
modified parton distribution functions (nPDFs) on x and 
squared momentum transfer, Q 2 , have been extracted 
by several groups from data that include deep inelas- 
tic electron- nucleus scattering (DIS) and Drell-Yan (DY) 
dilepton production in p + A collisions. The DIS and 
DY data together provide strong constraints on valence 
and sea quark modifications [2-5 . Including neutrino- 
induced DIS data from heavy targets can discriminate 
between quarks and antiquarks [3j 0] . Inclusive pion pro- 
duction data from RHIC have also been incorporated to 
better constrain the gluon modifications [4j [5] . 

The measurements used to extract the nPDFs cited 
above [2H5] were all averaged over impact parameter. 
Therefore these nPDFs represent the parton modifica- 
tion averaged over the entire nucleus. If the modification 
of these nPDFs is desired as a function of the impact pa- 
rameter, a specific dependence has been assumed [6]. A 
different approach, employing Gribov theory and incor- 
porating diffractive data, allows the spatial information 
to be retained [7]. 

Recently, the impact parameter dependence of the 
EPS09 [5 and EKS98 [8 nPDF's has been parameter- 
ized [9 using the target mass dependence of the EPS09 
and EKS98 parameter sets themselves. Terms up to 
fourth order in the nuclear thickness were necessary to 
produce A-independent coefficients. 



In this paper, we address the impact parameter de- 
pendence of gluon shadowing in a different way, using 
the collision centrality and rapidity dependence of the 
J yields measured in y/s NN = 200 GeV d+Au colli- 
sions at RHIC [10]. We were motivated by the observa- 
tion [lOj [11] that the onset of J / ijj suppression at forward 
rapidity suggests a quadratic or higher dependence on 
the nuclear thickness function at impact parameter tt, 
T A (r T ). 

Gluon-gluon interactions dominate J/ip production in 
high-energy hadronic collisions. Therefore, J/ip produc- 
tion in p(d)+A collisions must reflect the gluon modi- 
fication in the nuclear target. However, the measured 
modifications of J/ip yields in p(d)+A collisions relative 
to p + p collisions are also sensitive to the breakup of 
bound cc pairs by collisions with nucleons as the pairs 
pass through the medium. This effect, as well as effects 
due to any processes aside from shadowing, are usually 
parameterized by an effective absorption cross section, 
cr a b s , fitted to the measured data (see e.g. Ref. [12 ). 
The main goal of this work was to determine whether 
the effects of shadowing could be separated from those 
embodied in cr a b s . Because the magnitude of the effect 
due to cr a bs depends exponentially on nuclear thickness 
for a constant <r a b s , such separation may be possible if 
shadowing has a stronger thickness dependence. 

We assumed that the shadowing modification, inte- 
grated over all tt, could be described by the EPS09 gluon 
modification [5] and fitted the tt dependence of shadow- 
ing and the magnitude of cr a b s to the data. A Glauber 
Monte Carlo calculation [13] allows the modification, cal- 
culated for individual nucleon-nucleus collisions, to be 
correctly averaged and integrated over centrality, rapid- 
ity and pt- It also accounts for the effects of trigger 
efficiency in peripheral events. The Glauber parameters 
used here are identical to those used by PHENIX when 
calculating the experimental centrality distributions [10] . 
For that reason, we found we could drop the uncertain- 
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ties from the measured data points that are associated 
with the mean number of binary collisions because they 
are common to both the data and the calculations. The 
Woods-Saxon nuclear density distribution has a radius of 
6.34 fm and a diffuseness of 0.54 fm. No nuclear modifi- 
cation is assumed for the deuteron. The baseline J/ip pr 
and rapidity distributions used in the calculation were 
the p + p distributions measured by PHENIX [14] . 

The values of the target momentum fraction, #2, and 
squared momentum transfer, Q 2 , were assumed to obey 
approximate 2 1 kinematics as functions of J/ip ra- 
pidity and transverse momentum: 
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where M is the J/xJ; mass. The 2—^1 kinematics are not 
strictly correct since a high pt J / ip requires production of 
an associated hard parton. However, Eq. differs from 
exact 2 1 kinematics since the pr of the J/ip is finite. 
This approximation is close to the inclusive J/ip kinemat- 
ics in the CEM calculation described in Ref. [15] , NLO in 
the total cross section. Thus modifications of the gluon 
distribution in the nucleus are similar in the CEM calcu- 
lation and the Glauber Monte Carlo using Eq. 0. The 
results also agree with those found using PYTHIA [TT] , 

We tested several assumptions of the vt dependence 
of shadowing. Each of the postulated behaviors had one 
or two parameters that were adjusted to the data, along 
with the magnitude of cr a b s . The Glauber Monte Carlo 
was employed to compare the modification of the J/ip 
yields to the PHENIX d+Au data, averaged over four 
centrality bins and twelve rapidity bins, then integrated 
over all pt- We fitted <r a b s and the shadowing parameters 
to the PHENIX data using a modified x 2 function that 
properly accounts for all of the experimental uncertain- 
ties. Systematic uncertainties are included by moving 
the data points through ±3a in each systematic uncer- 
tainty while taking an appropriate x 2 penalty [16] . Indi- 
vidual parameter uncertainties were evaluated by finding 
the values at which x 2 increased by 1.0 with all other 
parameters re-optimized. 

We first assume that the shadowing strength is pro- 
portional to T%(tt), 



M shad = l-(l-R g (x,Q 2 )) 



a(n) 



(2) 



Here R g {x, Q 2 ) is the EPS09 gluon modification and the 
normalization factor, a(n), is adjusted so that the in- 
tegral over all impact parameters returns the average 
EPS09 modification. The power n was allowed to be 
unphysically large, n < 50, while constraining the mod- 
ification to be > 0. By allowing such arbitrarily large 
values of n, we test the sensitivity of the data to the 
centrality dependence of the shadowing. 



We began by fitting the values of both a a b s and n com- 
pletely independently at each rapidity. This is not quite 
realistic, because it ignores the correlations among some 
of the systematic uncertainties within each of the three 
spectrometer arms. However it provides an indication of 
how the centrality dependence varies with rapidity. The 
X 2 contours in a a b s and n corresponding to Ax 2 = 1.0 
and 2.3 are shown in Fig. [I] for the most backward ra- 
pidity, midrapidity, and the most forward rapidity. At 
midrapidity, the fits are insensitive to n because the shad- 
owing effects are weak. At the most forward and back- 
ward rapidities the optimum n is large, n > 10, indi- 
cating that the data require a strongly nonlinear onset 
of shadowing or antishadowing as a function of impact 
parameter. Additionally, there is relatively little correla- 
tion, and thus little ambiguity, between n and <r a b s for n 
greater than a few. 

The optimum value of a a b s and the corresponding un- 
certainty at each rapidity is shown by the red squares 
in Fig. [2] The <r a b s values are reasonably well defined, 
with a minimum near midrapidity. Because the data at 




FIG. 1: The x 2 distributions for the most backward, mid 
and most forward rapidities when cr a bs and n are optimized 
separately at each rapidity. 

all rapidities are consistent with n > 5, we repeated the 
fit assuming a rapidity-independent value of n. For this 
global fit, as for the one discussed later, those systematic 
uncertainties that are correlated within each experimen- 
tal arm were forced to move together. The optimum 
global power was n = 15^4- The x 2 /dof was 1.94 over- 
all. The best fit values of a a b s with n = 15 are given by 
stars in Fig. [2] 

The fits favor (or at midrapidity are consistent with) 
a shadowing modification that is negligible at large r^, 
but turns on sharply as vt decreases below ~ 2 — 3 fm. 
This behavior suggested that a step function onset of 
shadowing, including a radius, and a diffuseness, 
parameter, 
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would be more appropriate. 

Again, we began by fitting the parameters R and d. 
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FIG. 2: The optimum values and uncertainties on cr a bs as 
a function of rapidity obtained from the four fits described 
in the text. For clarity, the rapidities for each fit are slightly 
offset. 

along with <7 a b s , individually at each rapidity. The opti- 
mum cr a bs values, triangles in Fig. [2j are very similar to 
those obtained from the earlier fits employing T^r^). 
The fits favor R values of about half the Au radius, 
R < 3.5 fm, and a small diffuseness parameter although 
they are relatively insensitive to the value of d. 

Finally, we fitted cr a b s to the data at each rapidity while 
fitting global values of R and d. The x 2 /dof was 1.96 over 
all values of y. The best fit a a b s values at each rapidity 
are shown by circles in Fig. [2] The x 2 contours in R and 
d are shown in Fig. [3| The optimum global parameter 
values are R = 2.4±£|? fm and d = 0.12^;^ fm, where 
the uncertainties are obtained from the maximum extent 
of the Ax 2 = 1.0 contour. The fit results and RdAu are 
compared as a function of tt in Fig. [4} The dashed curves 
indicate the uncertainty in RdAu due to the uncertainty 
in cr a bs. Because x 2 includes the global uncertainties on 
the data, the best fit values may have a slight vertical 
offset from the data to achieve the best overall x 2 - 

The cr a b s values obtained from our various assumptions 
of the centrality dependence of shadowing are all com- 
pared in Fig. [2] The fitted values of cr a bs are well de- 
fined and show little dependence on the shadowing pre- 
scription, suggesting that the effects of shadowing can be 
separated from the effective breakup of the cc pair due 
to collisions with nucleons. Evidently, the strongly non- 
linear onset of shadowing with decreasing tt exhibits a 
distinctly different centrality dependence than the expo- 
nential dependence of <r a bs on ^a(Jt)- 

The gluon modification obtained from the EPS09 pa- 
rameterization using Eq. (|3|, with global fit values of R 
and d, is shown by the solid red line in Fig. [5] The effect 
of the combined uncertainty in R and d can be visual- 
ized by plotting the modifications for all combinations of 
R and d that produce a x 2 value inside the Ax 2 = 2.3 
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FIG. 3: The x 2 contours obtained from the fits with global 
values of R and d, Eq. (|3|. The uncertainties in R and d are 
taken from the maximum extent of the Ax 2 = 1 contour. 



contour [TT] (see Fig. J3|. These are represented by the 
thin blue lines in Fig. [5] In all cases, the calculated mod- 
ification is significant only for tt ^ 3 fm. Therefore, we 
conclude that the data place strong constraints on the 
nuclear modification so that it becomes significant only 
at small vt- The modification obtained with the best fit 
global power, T2(tt), n = 15, is shown as the solid or- 
ange line in Fig. [5] Although there is some difference in 
the details at small tt, albeit within the uncertainties, 
the two prescriptions give essentially the same value of 
X 2 /dof. Thus the data appear to be insensitive to the 
detailed shape of the modification at low ry, because the 
d+Au centrality bins are wide with significant overlap. 

We compare our results obtained from the fits to the 
J/ip data with those given by the newly- available impact 
parameter dependent EPS09s set [9]. It is shown by the 
dotted magenta line in Fig. |5j The EPS09s result has 
a much weaker dependence on tt than is obtained from 
our fits. A study of the target mass systematics of J/i/j 
production in d+^4 collisions at RHIC may shed light on 
the source of this pronounced difference. 

In summary, we have fitted the centrality and rapidity 
dependent PHENIX ^/s^ = 200 GeV d+Au J /if) data 
with Glauber calculations employing an effective breakup 
cross section, a a b s , with several prescriptions for the im- 
pact parameter dependence of the central EPS09 gluon 
shadowing parameterization. The fits properly account 
for all of the experimental systematic uncertainties. We 
find little ambiguity between a a b s and the functional form 
of the centrality dependence of shadowing. The values of 
<j a bs exhibit a characteristic rapidity dependence, with 
a minimum at midrapidity. The centrality dependence 
of shadowing suggested by the data turns on sharply for 
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FIG. 4: Comparison to the data of the best fits with global 
R and d values. 
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FIG. 5: The gluon modification from the best fit global R and 
d parameters (solid red line), along with the modifications 
from all combinations of R and d that fall within the A% 2 = 
2.3 contour in Fig. [3] (thin blue lines). The modification from 
the best fit global analysis of T2(tt) [n = 15) is shown by 
the solid orange line. The dashed magenta line is the recently 
released EPS09s impact parameter dependence [9]. 
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tt < 3 fm, in significant disagreement with the weaker tt 
dependence of EPS09s. Indeed, the EPS09s dependence 
is somewhat weaker than the linear dependence on the 
thickness function assumed in Ref. [6]. While we have 
employed only the central EPS09 set in our calculations, 
using all 31 EPS09 sets would not affect our overall con- 
clusions regarding the sharp turn on of shadowing with 
tt, only increase the uncertainty in the value of cr a b s as 
a function of rapidity. Rather, the strong impact param- 
eter dependence suggested here is in better accord with 
the 'hot spots' conjectured in a saturated medium of high 
gluon density. Such behavior at backward rapidity, in the 
antishadowing region, is, however, at odds with the sat- 
uration picture. 
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